/* Simple benchmark file, runs only one iteration for given size. Doesn't do
 * any of the fancy error checking stuff. Makes it easier to diagnose
 * performance stuff using gdb/other tools
 * It is highly recommended that you check the correctness of the algorithm
 * using the provided benchmark.c
 */


#define TEST_SIZE 768

#include <stdlib.h> // For: exit, drand48, malloc, free, NULL, EXIT_FAILURE
#include <stdio.h>  // For: perror
#include <string.h> // For: memset

#include <float.h>  // For: DBL_EPSILON
#include <math.h>   // For: fabs

#ifdef GETTIMEOFDAY
#include <sys/time.h> // For struct timeval, gettimeofday
#else
#include <time.h> // For struct timespec, clock_gettime, CLOCK_MONOTONIC
#endif

/* Your function must have the following signature: */
extern const char* dgemm_desc;
extern void square_dgemm (int, double*, double*, double*);

double wall_time ()
{
#ifdef GETTIMEOFDAY
  struct timeval t;
  gettimeofday (&t, NULL);
  return 1.*t.tv_sec + 1.e-6*t.tv_usec;
#else
  struct timespec t;
  clock_gettime (CLOCK_MONOTONIC, &t);
  return 1.*t.tv_sec + 1.e-9*t.tv_nsec;
#endif
}

void die (const char* message)
{
  perror (message);
  exit (EXIT_FAILURE);
}

void fill (double* p, int n)
{
  for (int i = 0; i < n; ++i)
    p[i] = 2 * drand48() - 1; // Uniformly distributed over [-1, 1]
}

/* The benchmarking program */
int main (int argc, char **argv)
{
  printf ("Description:\t%s\n\n", dgemm_desc);
  int test_size = TEST_SIZE; 

  /* allocate memory for all problems */
  double* buf = NULL;
  buf = (double*) malloc (3 * test_size * test_size * sizeof(double));
  if (buf == NULL) die ("failed to allocate largest problem size");

  /* Create and fill 3 random matrices A,B,C*/
  double* A = buf + 0;
  double* B = A + test_size*test_size;
  double* C = B + test_size*test_size;

  fill (A, test_size*test_size);
  fill (B, test_size*test_size);
  fill (C, test_size*test_size);

  /* Measure performance (in Gflops/s). */

  double Gflops_s, seconds = -1.0;

  seconds = -wall_time();
  square_dgemm (test_size, A, B, C);
  seconds += wall_time();

  /*  compute Mflop/s rate */
  Gflops_s = 2.e-9 * test_size * test_size * test_size / seconds;

  printf ("Size: %d\tGflop/s: %.3g\n", test_size, Gflops_s);

  free (buf);

  return 0;
}
